Noninvasive measurement of dissipation in colloidal systems 
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According to Harada and Sasa [Phys. Rev. Lett. 95, 130602 (2005)], heat production generated 
in a non-equilibrium steady state can be inferred from measuring response and correlation functions. 
In many colloidal systems, however, it is a nontrivial task to determine response functions, whereas 
details about spatial steady state trajectories are easily accessible. Using a simple conditional 
averaging procedure, we show how this fact can be exploited to reliably evaluate average heat 
production. We test this method using Brownian dynamics simulations, and apply it to experimental 
data of an interacting driven colloidal system. 
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Introduction. — The phenomenon of dissipation distin- 
guishes equilibrium from non-equilibrium systems. For 

macroscopic systems, dissipation can be inferred cither 
directly through measuring temperature changes or from 
the known external work applied to a system. At least for 
non-equilibrium steady states (NESSs) the latter is equal 
to the dissipation. For small systems such as colloidal 
particles or molecular motors, measuring dissipation is 
highly nontrivial. Calorimetric methods on the single 
particle or molecule level fail due to the tiny values of the 
heat generated by single degrees of freedom. In principle, 
the framework of stochastic thermodynamics [1-3] allows 
to apply the first law to phenomena on this scale from 
which the exchanged heat could be extracted if both, the 
externally applied work to drive the system and the inter- 
nal energy change could be measured, where the latter 
vanishes in a NESS. In practice, however, knowing the 
external force applied to a colloidal particle by a laser- 
field beyond the paradigmatic harmonic trap [4] is quite a 
challenge [5, 6]. Likewise, measuring directly the amount 
of adenosine triphosphate (ATP) molecules hydrolyzcd 
by a single moleciilar motor is impossible. Thus, the in- 
gredients of using the first law to infer dissipation, in 
general, arc not directly accessible. 

Harada and Sasa suggested an exact relation quan- 
tifying heat production in terms of the violation of the 
fluctuation-dissipation relation (FDR) for systems with a 
single degree of freedom obeying overdamped Langevin 
dynamics [7]. Extensions to many-body systems both 
for over- and imderdamped Langevin dynamics [8], as 
well as Hamiltonian systems [9] followed shortly after. 
Their results have been experimentally tested for an op- 
tically driven colloidal system [10] and applied to gain 
information about the non-equilibrium energetics of Fi- 
ATPase [11]. Further generalizations and experiments 
applying this technique comprise Langevin systems in- 
cluding memory [12, 13], quantum Langevin dynam- 
ics [14], and FDR violations involving field variables [15]. 
For molecular motors connected to a colloidal probe, an- 
other approach has been put forward, from which infor- 



mation about dissipation can be deduced [16, 17]. 

Specifically, the Harada-Sasa relation expresses the av- 
erage heat production rate as 



Qi(u;)-2Ti?^(u;) 



dtu 
2^ 



(1) 

where Boltzmann's constant is set to unity, T is the tem- 
perature of the surrounding heat bath and /i,; is the mo- 
bility of the ith of the n degrees of freedom. The response 
function 



Rij {t) = 



(2) 



quantifies the change of the velocity Xi to a force per- 
turbation applied to the jth degree of freedom. Further- 



more, 



Qj{t)^{{xi{t)-v^,)ixj{0)-v'.)) 



(3) 



is the velocity autocorrelation function, and g{u)) = 
g(t) exp [iu}t)dt / {27^) is the Fourier transform of an 
arbitrary function g{t). The prime in Eq. (1) denotes the 
real part of a complex-valued function. 

This exact relation yields insight into the origins of 
FDR violations and can be used readily if response func- 
tions are accessible. Determining response fimctions for 
driven colloidal systems, however, often poses a nontriv- 
ial and time-consuming task due to the need to perturb 
the system from its steady state. The perturbation must 
be sufficiently small to stay within the linear response 
regime and has to be applied to each degree of freedom 
separately. Determining response functions in frequency 
space requires separate measurements for each frequency, 
which must be done up to sufficiently high frequencies to 
make sure that the integral in Eq. (1) converges. If deter- 
mined in temporal space, the system must be observed 
during its relaxation, which excludes the possibility to 
average over the time coordinate, thus substantially in- 
creasing the statistical effort. 
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In experiments and simulations, good statistics for tra- 
jectories in a NESS is readily available. In the follow- 
ing, we will employ a conditional averaging procedure on 
these trajectories from which the mean local velocity field 
can be obtained. Combined with the measured station- 
ary distribution, this field then yields the average heat 
production rate. This procedure avoids difficulties aris- 
ing from the application of external perturbations and 
reduces the effort to create the statistics needed. It thus 
complements the approach given by Harada and Sasa by 
making the average heat production rate easily accessi- 
ble from steady state trajectories. Likewise, it does not 
require any knowledge about the applied forces or inter- 
actions within the system. 

Theory. — We introduce this method for systems whose 
dynamics is governed by a set of coupled overdamped 
Langevin equations 

X = AiF(x) + C, (4) 

with the coordinate x = (xi, . . . ,a;„), and the force F 
composed of conservative and nonconservative contribu- 
tions. Interactions with the surrounding solvent are mod- 
eled by Gaussian white noise ^ with zero mean (C(^)) = 
and correlations 

mC'^{t'))=2fiTd{t-r). (5) 

The stationary probability distribution function p^(x) 
follows from the Smoluchowski equation according to 

= -V-j^(x) (6) 

with the probability current 

(x) = fiF (x)p^ (x) - Vp^ (x) . (7) 

The closely related mean local velocity is defined as the 
conditional average of the fluctuating velocity x at posi- 
tion X, 

u%x) = lim (x(i + At)-x(t-At)|x(t) = x)/(2At). 

(8) 

Here, the Stratonovich convention has to be employed, 
thus evaluating the spatial variable in mid-step position. 
For the remainder of this Rapid Communication, we will 
always assume this convention. This relation enables us 
to evaluate the mean local velocity from NESS trajec- 
tories. Carrying out the average analytically, one ob- 
tains [1] 

iv^(x) = j^(x)//(x) = /iF(x) - ^TVlnp^(x). (9) 

The average heat production rate along a stochastic 
trajectory x{t) is given by [2] 

(g>^(F(x(i))-i(t)). (10) 



Since x is a stochastically fluctuating quantity some care 
needs to be taken while averaging. In the Stratonovich 
scheme, one has [1] 

(g)-(F(x).z.^(x)). (11) 

Now, all fluctuating quantities have been replaced. Using 
Eq. (9) we substitute F and obtain 

(q) = /^-i J dxi/"(x)V(x) +T J dxiy"(x) • Vp"(x). 

(12) 

The second integral on the right-hand side vanishes, since 

J dxv^{x) ■ Vp"(x) = J fixj"(x) • Vlnp"(x) 

= ^ydx[V-j^(x)]lnp^(x) = 0. (13) 

Here, we have used Eq. (6) after a partial integration. 
The boundary term vanishes due to the periodicity of 
j*' and p^. Therefore, the average heat production rate 
reads 

{q) = f^-' J dxi.^(x)V(x) = ^r'{l.^{xf). (14) 

Given a long trajectory x(t), we can determine the sta- 
tionary distribution and via Eq. (8) also the mean 
local velocity i/^. Therefore, we are able to evaluate the 
average heat production rate in the system simply by 
recording particle trajectories in the NESS. 

In the following, we will first illustrate the validity and 
usefulness of the method by using Brownian dynamics 
simulations to model a driven system consisting of two 
coupled colloidal particles. Second, we apply the method 
to data for the experimental system and determine the 
average heat production rate as a function of the coupling 
strength. 




FIG. 1: Schematic representation of the system: Two para- 
magnetic colloidal particles driven along two rings of radius 
R by constant forces fi. The position Xi of the ith particle 
is the arc length measured in the counterclockwise direction. 
The magnetic field applied perpendicular to the rings induces 
a magnetic moment in each particle. The strength of the re- 
sulting repulsive interaction is quantified by the dimensionless 
plasma parameter F. 
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FIG. 2: Comparison between the force fields F(x) for the 
two-ring system obtained (a) analytically and (b) from NESS 
trajectories via Eq. (9) for parameter set I [see Table I]. 



Simulation. — We use Brownian dynamics simulations 
to test this method for the following, experimentally ac- 
cessible [18], colloidal system. Two paramagnetic col- 
loidal particles are driven on two nonoverlapping rings 
of radius R [see Fig. 1]. By applying a small homoge- 
neous magnetic field B normal to the plane containing 
the rings, the particles acquire parallel magnetic dipole 
moments m w aB with a ~ 5.9 x 10~^^ Am^/T. The re- 
sulting repulsive interaction is described by the potential 
T4^(x) — ^o"T.^/(47rr'^(x)), where is the magnetic con- 
stant and r(x) the distance between the particles. In or- 
der to quantify the interaction strength, we introduce the 
dimensionless plasma parameter F = AVF/T, where /SW 
is the difference between the maximum and minimum in 
the interaction energy. Apart from the interaction, the 
total force acting on the particle on the ith ring, 



F,(x) = /,-a,^ {y^{x,) + w{^)), 



(15) 



contains the sinusoidal potential Vi = Ai am.{xi/ R—4>i) of 
amplitude Ai and phase shift (pi , and the constant driving 
force fi. Using these forces, the dynamics of the system 
is given by the Langevin equation (4). The equations of 
motion are integrated via a stochastic Runge-Kutta algo- 
rithm [19] using a time step of At = 0.0001 (^m)V(AiT). 

The most crucial point of our method is the determi- 
nation of the mean local velocity ly^ and the probability 
distribution function p^. In order to check this point, we 
proceed as follows. First, we compute the force field F 
analytically, which is possible since in this simulation all 
potential parameters are known in detail. Second, we 
compute F via Eq. (9) by evaluating both and the 
gradient field of from the simulated data. In Fig. 2, 
we compare the differently obtained force fields and find 
very good agreement. 

We calculate the average heat production rate using 
two different methods: (A) We use the trajectory-based 
method introduced in this Rapid Communication, i. e.,we 
determine v'^ and from the simulated trajectories and 
use Eq. (14) to obtain (g). (B) As a test, we determine 
(F • x) along the complete trajectory for each particle. 
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FIG. 3: (a): Average heat production rates (g) determined 
from methods (A) and (B) as described in the text for three 
difi'erent parameter sets I, II, and III [ see Table I]. To relate 
(q) to an experimental scale we choose a — 5.2 \Lm. Error bars 
are smaller than 1%. (b)-(d): Corresponding evaluated mean 
local velocity fields i-'°(x). Regions without arrows are rarely 
visited by the system. 
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h [kBT/tmr] 
/2 [keT/^m] 
v4i [ksT] 

A2 [kBTl 



R [\im\ 



310.0 
50.0 
-65.0 
225.0 
175.0 

-7r/5 
-37r/5 
3.5 



670.0 
56.0 
-54.0 
175.0 
170.0 
-7r/5 
-37r/5 
3.5 



310.0 
25.0 
-35.0 
175.0 
175.0 
-7r/5 
-37r/5 
3.5 



center-to-center distance [yan\ 17.0 17.0 17.0 



TABLE I: Parameter sets /, //, III. 



The sum then yields the true average heat production 
rate in this simulation run. Of course, this procedure is 
possible only if details of the interaction are known, as is 
the case in these Brownian dynamics simulations. 

In Fig. 3, we present a comparison between the differ- 
ently evaluated average heat production rates determined 
for three different parameter sets [see Table I]. 

We find very good agreement, supporting the consis- 
tency between the present approach (A) and the reference 
result (B). The most critical point in the application of 
method (A) is the determination of in those phase- 
space regions where the system is found with low proba- 
bility. The average heat production rate (g), however, is 
computed from an integral [see Eq. (14)] in which {v^Y 
enters weighted by p^ . Since in rarely visited regions this 
weight is small, their contributions to {q) should play 
a minor role, provided that stays finite. This is the 
case in our system, and should hold for many other well- 
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behaved systems as well. Therefore, to obtain averaged 
quantities such as (q) , only regions in configuration space 
of substantial weight must be sampled accurately. Re- 
gions of minor probability will only add small-sized cor- 
rections to the result. 

Experiment. — Having demonstrated its validity, we 
now apply the method to the experimental system [18]. 
Two paramagnetic colloidal particles with diameters a — 
5.2 |j.m are driven by a scanning laser beam [20, 21] along 
two rings with radii R = 3.5 |i.m and a center-to-center 
distance of 17 |J.m. We checked for hydrodynamic interac- 
tions by measuring the mean velocity of a particle driven 
along a circle with constant potential in the presence of 
another particle at minimal distance d. For d ^ 9 \xni, we 
do not find any deviation from the d — >■ oo limit. We thus 
conclude that hydrodynamic interactions are negligible in 
this case [23]. Using digital video microscopy, we track 
the NESS trajectories xi{t) and X2{t) with a spatial and 
temporal resolution of 20 nm and 25 ms, respectively. 
We determine the average heat production rate for dif- 
ferent magnetic fields, B ^ 40 mT, thus changing the 
plasma parameter in the range of ^ F ^ 1100. 

In the absence of coupling, the particles move inde- 
pendently along their tilted potentials with a mean cir- 
culation time of 12 s. This motion leads on average to a 
heat production rate of (q) = fj,~^ {v^ (x.)'^) = 199 keT/s 
with = Snrja the bare mobility of one particle and rj 
the viscosity of the solvent. Under strong coupling con- 
ditions, (q) reaches higher values [see Fig. 4 (a)]. This 
fact becomes obvious when focusing on the correspond- 
ing mean local velocity fields i/^. For F = 0, i/*^ shows 
a large region with very small velocities [see Fig. 4 (b)] 
corresponding to potential minima at Xi/R = 7r/4 and 
X2/ R = 2>-K /4:. For F ^ 480, the motion of the two parti- 
cles synchronizes due to the repulsive interaction as fol- 
lows [see Fig. 4 (c)]. While the first particle fluctuates 
around its potential minimum, the second one completes 
one circulation. Arriving at the point where the parti- 
cles are closest, this particle pushes the first one over 
its barrier. Then, the roles switch and the first particle 
completes a circulation while the second one fluctuates 
around its potential minimum. This effect shortens the 
time to overcome the potential barrier considerably for 
both particles and results, on average, in a faster circula- 
tion, and thus a larger (g) than in the uncoupled case. In 
the intermediate regime, the particles may still overtake 
one another. During the overtaking process, the faster 
particle hinders the other one in overcoming the potential 
barrier, effectively slowing down its motion. This on av- 
erage slower motion produces less heat, and consequently 
(g) is non-monotonic in F [see Fig. 4 (a)]. We emphasize 
that we obtained robust results for {q) for NESS trajec- 
tories with a length of only 30 min corresponding to 150 
full revolutions of each particle. 

Concluding perspectives. — In this Rapid Communica- 
tion, we have presented an approach to determine the 
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FIG. 4: (a): Average heat production rates (g) for different 
plasma parameters F. The parameters for the NESSs are /i = 
56 ksT/^m, Ai = 175 keT, and /2 = -51 keT/nm, Ai = 
164 keT obtained by the method presented in Ref. [22]. Error 
bars are obtained by evaluating the data set for four parts of 
equal size and computing the standard deviation, (b)-(c); 
Corresponding mean local velocity fields i''^(x) determined 
from the NESS trajectories for plasma parameters (b) F = 0, 
and (c) r = 1150, respectively. 



average heat production rate for colloidal systems in a 
NESS. While the Harada-Sasa method quantifies dissipa- 
tion via response and correlation functions, the present 
method solely uses information obtained from steady 
state trajectories without the need to perturb the system. 
Neither does the method require any information about 
the underlying potentials or driving forces. Recording 
particle trajectories in the NESS is sufficient. This ap- 
proach thus constitutes a complementary way to measure 
dissipation, which is easily implemented particularly in 
set-ups concerned with small colloidal systems. By re- 
placing the stationary distribution function and the sta- 
tionary mean local velocity field by their time-dependent 
counterparts [1], a generalization of the present method 
should also be applicable to relaxing and time-dependent 
systems. 
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